Jednowymiarowa metoda elementów skończonych
Dla danej siatki obliczeniowej \( \{ T_{hp }= \left(K, X\left(K\right), \Pi_p \right) \}_K \) siatka referencyjna to rodzina \( T_{\frac{h}{2}p+1} \) elementów skończonych
\( T_{\frac{h}{2}p+1}=\{ \left( K, X\left(K\right), \Pi_p \right) \}_K { } \) takich że
\( \forall K \in T_{ \frac{h}{2}p+1} \exists K_1,K_2 \in {\cal P}(T_{hp },K) \) takie że
\( K=K_1 \cup K_2, \textrm{meas}K_1\cap K_2=0 \),
\( X(K_1),X(K_2) \in {\cal P}(T_{hp },X(K)), dim X(K)=dimX(K_1)+1=dimX(K_2)+1 \) gdzie \( {\cal P } (T_{hp},K ) \) i \( {\cal P }(T_{hp},X(K)) \) oznaczają projekcje na pierwszy i drugi komponent \( \left( K, X\left( K \right), \Pi_p \right) \). Wzór wymiaru przestrzeni funkcji kształtu na nowych dwóch elementach wynika z dodania jednej nowej funkcji kształtu, stosownie poprzez zwiększenie stopnia wielomianu o 1 we wnętrzu elementu.
Przestrzeń aproksymacyjna nad jednowymiarową siatką obliczeniową \( T_{hp}=\{ \left(K, X\left(K\right), \Pi_p \right) \}_K \) to
\( V_{hp }= \{ v \in C(\Omega): \forall K \in {\cal P }( T_{hp},K): {\cal P }( v,K) \in X(K) \} \)
gdzie \( {\cal P }( T_{hp},K) \) to zbiór przedziałów reprezentujących geometrię elementów wyciągnięty z trójki reprezentującej jednowymiarową siatkę obliczeniową, \( {\cal P }( v,K) \) to projekcja funkcji na przedział reprezentujący geometrię elementu.
Definicja 4: Odwzorowania definiujące składanie lokalnych funkcji kształtu na globalne funkcje bazowe nad siatką obliczeniową
Niech \( T_{hp}=\{ \left(K, X\left(K\right), \Pi_p \right) \}_K \) oznacza jednowymiarową siatkę obliczeniową.
Niech \( \{ e_i^{hp} \}_i \) oznacza bazę przestrzeni \( V_{hp }= span\{ e_i^{hp} \} \).
Niech \( \chi^K_k \in X \left( K\right) \) oznacza funkcje kształtu nad elementem \( K \).
Wówczas \( \forall K \in {\cal P }( T_{hp},K), \forall i, \exists k : {\cal P}(e_i^{hp},K) = \chi^K_k \).
Istnieje odwzorowanie odwrotne \( {\cal I}^2 \ni (k,K)\rightarrow i(k,K)\in {\cal I} \) które przypisuje numer
\( i(k,K) \) globalnej i-tej funkcji bazowej związanej z lokalną \( k \)-tą funkcją kształtu nad elementem \( K \).
W przypadkach rozpatrywanych w tym podręczniku, odwzorowanie to jest izomorficzne.
Definicja 5: Przestrzeń aproksymacyjna nad jednowymiarową siatką referencyjną
Przestrzeń aproksymacyjna nad jednowymiarową siatką obliczeniową \( T_{\frac{h}{2}p+1}=\{ \left(K, X\left(K\right), \Pi_p \right) \}_K \) to
\( V_{\frac{h}{2}p+1 }= \{ v \in C(\Omega): \forall K \in {\cal P }( T_{\frac{h}{2}p+1},K): {\cal P }( v,K) \in X(K) \} \)
gdzie \( {\cal P }( T_{\frac{h}{2}p+1},K) \) to zbiór przedziałów reprezentujących geometrię elementów wyciągnięty z trójki reprezentującej jednowymiarową siatkę obliczeniową, \( {\cal P }( v,K) \) to projekcja funkcji na przedział reprezentujący geometrię elementu.
Definicja 6: Odwzorowania definiujące składanie lokalnych funkcji kształtu na globalne funkcje bazowe nad siatką referencyjną
Niech \( T_{\frac{h}{2}p+1}=\{ \left(K, X\left(K\right), \Pi_p \right) \}_K \) oznacza jednowymiarową siatkę referencyjną.
Niech \( \{ e_i^{\frac{h}{2}p+1} \}_i \) oznacza bazę przestrzeni \( V_{\frac{h}{2}p+1 }= span\{ e_i^{\frac{h}{2}p+1} \} \).
Niech \( \chi^K_k \in X \left( K\right) \) oznacza funkcje kształtu nad elementem \( K \).
Wówczas \( \forall K \in {\cal P }( T_{\frac{h}{2}p+1},K), \forall i, \exists k : {\cal P}(e_i^{\frac{h}{2}p+1},K) = \chi^K_k \).
Istnieje odwzorowanie odwrotne \( {\cal I}^2 \ni (k,K)\rightarrow i(k,K)\in {\cal I} \) które przypisuje numer
\( i(k,K) \) globalnej i-tej funkcji bazowej związanej z lokalną \( k \)-tą funkcją kształtu nad elementem \( K \).
W przypadkach rozpatrywanych w tym podręczniku, odwzorowanie to jest izomorficzne.
Siatka referencyjna używana jest do szacownia błędu względnego na siatce obliczeniowej. Często kiedy opisuje się siatkę obliczeniową w kontekście siatki referencyjnej, nazywa się ją siatką rzadka, a siatkę referencyjną siatką gęstą, powstaje ona bowiem poprzez połamanie elementów na mniejsze elementy i zwiększenie stopnia aproksymacji wielomianowej o jeden.
Definicja 7: Problem metody elementów skończonych nad jednowymiarową siatką obliczeniową
Dla danej siatki obliczeniowej \( T_{hp}=\{ \left(K, X\left(K\right), \Pi_p \right) \}_K \) znaleźć współczynniki \( \{ u^{hp}_i \}_{i=1,...,N^{hp}} \) rozwiązania przybliżonego \( V \supset V_{hp } \ni u_{hp }=\sum_{i=1,...,N^{hp } } u_i^{hp } e_i^{hp } \) takie że \( \sum_{i=m,...,N^{hp } } u_m^{hp } B(e_m^{hp},e_n^{hp})=L(e_n^{hp}), n=1,...,N^{hp} \) gdzie \( B(e_m^{hp},e_n^{hp})=\int_0^l \left( a(x) \frac{de_m^{hp}(x)}{dx} \frac{de_n^{hp}(x)}{dx} +b(x)\frac{de_m^{hp}(x)}{dx}e_n^{hp}(x)+c(x)e_m^{hp}(x)e_n^{hp}(x)\right)dx + \beta e_m^{hp}(l)e_n^{hp}(l) \) oraz \( L(e_n^{hp})=\int_0^l f(x) e_n^{hp}(x)dx + \gamma e_n^{hp}(l) \). Baza \( \{ e^{hp}_i \}_{i=1,...,N^{hp}} \) przestrzeni aproksymacyjnej \( V_{hp} \) uzyskana jest poprzez sumowanie funkcji kształtu z przestrzeni \( X\left(K\right) = {\cal P}(T_{hp},X(K)) \) dla poszczególnych elemetów z siatki \( T_{hp} \) w globalne funkcje bazowe.
Definicja 8: Problem metody elementów skończonych nad jednowymiarową siatką referencyjną
Dla danej siatki obliczeniowej \( T_{\frac{h}{2}p+1}=\{ \left(K, X\left(K\right), \Pi_p \right) \}_K \) znaleźć współczynniki \( \{ u^{\frac{h}{2}p+1}_i \}_{i=1,...,n^{\frac{h}{2}p+1}} \) rozwiązania przybliżonego \( V \supset V_{\frac{h}{2}p+1 } \ni u_{\frac{h}{2}p+1 }=\sum_{i=1,...,N^{\frac{h}{2}p+1 } } u_i^{\frac{h}{2}p+1 } e_i^{\frac{h}{2}p+1 } \) takie że \( \sum_{m=1,...,N^{\frac{h}{2}p+1 } } u_m^{\frac{h}{2}p+1 } B(e_m^{\frac{h}{2}p+1},e_n^{\frac{h}{2}p+1})=L(e_n^{\frac{h}{2}p+1}), n=1,...,N^{hp} \) gdzie \( B(e_m^{\frac{h}{2}p+1},e_n^{\frac{h}{2}p+1})=\int_0^l \left( a(x) \frac{de_m^{\frac{h}{2}p+1}(x)}{dx} \frac{de_n^{\frac{h}{2}p+1}(x)}{dx} +b(x)\frac{de_m^{\frac{h}{2}p+1}(x)}{dx}e_n^{\frac{h}{2}p+1}(x)+c(x)e_m^{\frac{h}{2}p+1}(x)e_n^{\frac{h}{2}p+1}(x)\right)dx +\\+ \beta e_m^{\frac{h}{2}p+1}(l)e_n^{hp}(l) \) oraz \( L(e_n^{\frac{h}{2}p+1})=\int_0^l f(x) e_n^{\frac{h}{2}p+1}(x)dx + \gamma e_n^{\frac{h}{2}p+1}(l) \).
Baza \( \{ e^{\frac{h}{2}p+1}_i \}_{i=1,...,N^{\frac{h}{2}p+1}} \) przestrzeni aproksymacyjnej \( V_{hp} \) uzyskana jest poprzez sumowanie funkcji kształtu z przestrzeni \( X\left(K\right) = {\cal P}(T_{\frac{h}{2}p+1},X(K)) \) dla poszczególnych elemetów z siatki \( T_{\frac{h}{2}p+1} \) w globalne funkcje bazowe.
Uwaga 1: Rozwiązywanie problemu metody elementów skończonych nad jednowymiarową siatką obliczeniową
Współczynniki \( \{ u^{hp}_i \}_{i=1,...,N^{hp}} \) oblicza się poprzez rozwiązanie układu równań \( {\bf A}{\bf u}={\bf b} \) gdzie
\( {\bf A}= \begin{bmatrix} B(e_1^{hp},e_1^{hp}) & \cdots & B(e_{N^{hp}}^{hp},e_1^{hp}) \\ \vdots & \cdots & \vdots \\ B(e_1^{hp},e_{N^{hp}}^{hp}) & \cdots & B(e_{N^{hp}}^{hp},e_{N^{hp}}^{hp}) \end{bmatrix} \).
\( {\bf u} = \begin{bmatrix} u_1^{hp} \\ \vdots\\ u_{N^{hp}}^{hp} \end{bmatrix} \).
\( {\bf b} = \begin{bmatrix} L(e_1^{hp}) \\ \vdots \\ L(e_{N^{hp}}^{hp}) \end{bmatrix} \).
Warunek brzegowy Dirichleta uzyskuję się poprzez wyzerowanie wierszy odpowiadającym funkcjom bazowym zdefiniowanym na brzegu Dirichleta, umieszczenie jedynki na przekątnej i wyzerowaniu prawej strony. Oddziaływanie warunku brzegowego na nasz układ równań zostało przeniesione na prawą stronę poprzez operację przesunięcia.
Całki \( B(e_m^{hp},e_n^{hp}) \) oraz \( L(e_n^{hp}) \) oblicza się element po elemencie stosując stosowne kwadratury Gaussa.
\( B(e_m^{hp},e_n^{hp})=\sum_k \left( a(x^k) \frac{de_m^{hp}(x^k)}{dx} \frac{de_n^{hp}(x^k)}{dx} +b(x^k)\frac{de_m^{hp}(x^k)}{dx}e_n^{hp}(x^k)+c(x^k)e_m^{hp}(x^k)e_n^{hp}(x^k)\right)Jac(x^k)w^k +\\+ \beta e_m^{hp}(l)e_n^{hp}(l) \) oraz \( L(e_n^{hp})=\sum_k f(x^k) e_n^{hp}(x^k) Jac(x^k)w^k + \gamma e_n^{hp}(l) \)
gdzie \( Jac(x^k) \) to wartość w punkcie kwadratury Gaussa Jakobianu odwzorowania z elementu wzorcowego na dany element.
Generację układu równań dla danej siatki obliczeniowej \( T_{hp}=\{ \left(K, X\left(K\right), \Pi_p \right) \}_K \) uzyskać można następującym algorytmem:
1 for \( K \in {\cal P}(T_{hp},K) \) (pętla po elementach siatki obliczeniowej)
2 for \( \psi_k \in X(K) \) (pętla po funkcjach kształtu na elemencie \( K \) )
3 \( {\bf b}(i(k,K)) += L(\psi_k) \) (obliczanie lokalnej kontrybucji całki na elementcie \( K \) )
4 for \( \psi_l \in X(K) \) (pętla po funkcjach kształtu na elemencie \( K \) )
5 \( {\bf A}(i(k,K),i(l,K)) += B(\psi_k,\psi_l) \) (obliczanie lokalnej kontrybucji całki na elementcie \( K \) )
Uwaga 2: Rozwiązywanie problemu metody elementów skończonych nad jednowymiarową siatką referencyjną
Współczynniki \( \{ u^{\frac{h}{2}p+1}_i \}_{i=1,...,n^{\frac{h}{2}p+1}} \) oblicza się poprzez rozwiązanie układu równań \( {\bf A }{\bf u} = {\bf b} \)
\( {\bf A} = \begin{bmatrix} B(e_1^{\frac{h}{2}p+1},e_1^{\frac{h}{2}p+1}) & \cdots & B(e_{N^{\frac{h}{2}p+1}}^{\frac{h}{2}p+1},e_1^{\frac{h}{2}p+1}) \\ \vdots & \cdots & \vdots \\ B(e_1^{\frac{h}{2}p+1},e_{N^{\frac{h}{2}p+1}}^{\frac{h}{2}p+1}) & \cdots & B(e_{N^{\frac{h}{2}p+1}}^{\frac{h}{2}p+1},e_{N^{\frac{h}{2}p+1}}^{\frac{h}{2}p+1}) \end{bmatrix} \).
\( {\bf u} = \begin{bmatrix} u_1^{\frac{h}{2}p+1} \\ \vdots\\ u_{N^{\frac{h}{2}p+1}}^{\frac{h}{2}p+1} \end{bmatrix} \).
\( {\bf b} =\begin{bmatrix} L(e_1^{\frac{h}{2}p+1}) \\ \vdots \\ L(e_{N^{\frac{h}{2}p+1}}^{\frac{h}{2}p+1}) \end{bmatrix} \).
Całki \( B(e_m^{\frac{h}{2}p+1},e_n^{\frac{h}{2}p+1}) \) oraz \( L(e_n^{hp}) \) oblicza się element po elemencie stosując stosowne kwadratury Gaussa.
\( B(e_m^{\frac{h}{2}p+1},e_n^{\frac{h}{2}p+1})=\sum_k \left( a(x^k) \frac{de_m^{\frac{h}{2}p+1}(x^k)}{dx} \frac{de_n^{\frac{h}{2}p+1}(x^k)}{dx} +b(x^k)\frac{de_m^{\frac{h}{2}p+1}(x^k)}{dx}e_n^{\frac{h}{2}p+1}(x^k)+c(x^k)e_m^{\frac{h}{2}p+1}(x^k)e_n^{\frac{h}{2}p+1}(x^k)\right) *\\ *Jac(x^k)w^k + \beta e_m^{\frac{h}{2}p+1}(l)e_n^{\frac{h}{2}p+1}(l) \) oraz \( L(e_n^{\frac{h}{2}p+1})=\sum_k f(x^k) e_n^{\frac{h}{2}p+1}(x^k) Jac(x^k)w^k + \gamma e_n^{\frac{h}{2}p+1}(l) \)
gdzie \( Jac(x^k) \) to wartość w punkcie kwadratury Gaussa Jakobianu odwzorowania z elementu wzorcowego na dany element.
Generację układu równań dla danej siatki obliczeniowej \( T_{\frac{h}{2}p+1}=\{ \left(K, X\left(K\right), \Pi_p \right) \}_K \)
uzyskać można następującym algorytmem:
1 for \( K \in {\cal P}(T_{\frac{h}{2}p+1 },K) \) (pętla po elementach siatki obliczeniowej)
2 for \( \psi_k \in X(K) \) (pętla po funkcjach kształtu na elemencie \( K \) )
3 \( {\bf b}(i(k,K)) += L(\psi_k) \) (obliczanie lokalnej kontrybucji całki na elementcie \( K \) )
4 for \( \psi_l \in X(K) \) (pętla po funkcjach kształtu na elemencie \( K \) )
5 \( {\bf A}(i(k,K),i(l,K)) += B(\psi_k,\psi_l) \) (obliczanie lokalnej kontrybucji całki na elementcie \( K \) )